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A  Candidate  Mesoscale 
Numerical  Cloud/Precipitation  Model 


1.  INTRODUCTION 

The  United  States  Air  Force  is  proceeding  with  a  program  called  the  Auto¬ 
mated  Weather  Distribution  System  (AWDS)  in  order  to  modernize  base -weather 
operations  of  the  Air  Weather  Service  through  the  use  of  advanced  computers  and 
peripheral  equipment.  This  modernization  generally  will  result  in  efficient  and 
accurate  analyses  and  forecasts  of  weather  systems.  In  particular,  the  primary 
objectives*  for  AWDS  are  to  provide:  (1)  operational  capability  for  better  detec¬ 
tion  and  transmission  of  critical  weather  elements;  (2)  modern  display  (e,  g.  , 
CRT)  and  dissemination  of  weather  observations;  (3)  a  modular  system  to  allow 
selective  combinations  of  components;  and  (4)  the  capability  for  expansion  and 
incorporation  of  new  display  and  sensing  techniques.  In  addition,  this  modern¬ 
ization  program  will  require  that  AWS  stations  have  the  capability  to  produce 
automated,  numerical  forecasts  of  relevant  weather  in  the  zero  to  12 -h  period 
and  for  spatial  domains  £400  km  (meso-beta  or  "theater"  scale).  The  motiva¬ 
tion  for  the  last  requirement  is  the  desire  to  employ  objective  numerical  weather 
prediction  techniques  that  produce  deterministic  mesoscale  forecasts  rather  than 


(Received  for  publication  9  October  198-5) 

1.  Air  Weather  Service  (1984)  C^abilities  Master  Plan  (CMP)  1984  -1985, 
Headquarters  AWS,  Scott  AFB,  Ill.  " 


relying  solely  on  subjective  inference  from  the  present  system  of  synoptic -scale 
guidance,  climatology,  terrain,  and  meteorological  observations.  However, 
adequate  numerical  treatment  of  the  dynamic  processes  that  occur  on  these  scales 
has  not  yet  been  realized. 

Improvement  in  mesoscale  forecast  skill  is  inhibited  by  problems  in  four 
principal  areas.  ^  The  first  area  concerns  numerical  and  practical  difficulties 
associated  with  realistic  treatment  of  lateral  and  lower  boundary  conditions  in 
the  model.  The  second  is  the  introduction  of  significant  error  in  the  analyzed 
initial  conditions  for  the  model  due  to  inadequate  spatial  and  temporal  resolution 
of  initial  upper-air  data  from  the  present  observational  network.  The  third  area 
is  the  difficulty  stemming  from  the  parameterization  in  the  model  of  the  sub -grid 
scale  physical  processes.  The  fourth  problem  area  is  the  Inherent  limit  in  the 
predictability  of  smaller -scale  circulations. 

The  Air  Force  Geophysics  Laboratory  (AFGL)  has  initiated  a  program  to 
develop  predictive,  sub-synoptic  scale  regional  models  specifically  designed  to 
provide  1-  to  2-h  forecasts  of  clouds  and  precipitation.  The  program  at  AFGL 
is  intended  to  select,  test,  and  evaluate  a  model  (or  models)  that  ultimately  could 
be  used  in  two  operational  modes:  (1)  by  nesting  within  a  larger-scale  model,  or 
(2)  on  a  stand-alone  basis  in  a  base-weather  station,  mini-computer  environment. 
The  model^  selected  for  these  studies  was  developed  at  the  National  Oceano¬ 
graphic  and  Atmospheric  Administration  (NOAA).  It  will  be  described  in  Section  2. 

While  much  of  the  development  effort  of  the  model  (both  nested  and  stand¬ 
alone  modes)  was  being  carried  out  by  contractors  at  the  NOAA  Environmental 
Research  Laboratory,  AFGL  scientists  concluded  a  baseline  test  for  a  prelim¬ 
inary  assessment  of  its  capabilities,  computer  requirements,  and  possible  weak¬ 
nesses.  This  technical  report  will  describe  the  NOAA  model,  the  problems  and 
weaknesses  discovered  during  the  baseline  test,  and  a  series  of  suggested  reme¬ 
dies  to  solve  these  problems. 


2.  Anthes,  R.  A.  (1983)  Regional  models  of  the  atmosphere  in  middle  latitudes, 

Mon.  Wea.  Rev.  Ill:  1306-  1335. 

3.  Nickerson,  E.  C.  (1979)  On  the  numerical  simulation  of  airflow  and  clouds 

over  mountainous  terrain,  Beit,  zur  Physik  der  Atmos.  52 : 16 1  -  177 . 

4.  Nickerson,  E.C.  ,  Richard,  E.,  Rosset,  R. ,  and  Smith,  D.  (1985)  On  the 

numerical  simulation  of  clouds,  rain,  and  airflow  over  the  Vosges  and 
Black  Forest  Mountains;  A  meso-|3  model  with  parameterized  micro¬ 
physics,  Mon.  Wea.  Rev,  (in  press). 


2.  DKSCKII'TION  OK  TIIK  MODKL 


2.1  Dyiiuinu' 

The  NOAA  model^  is  a  3D,  hydrostatic,  primitive  equation  model  based,  in 

part,  on  the  Anthes  and  Warner’  model.  Finite  difference  approximations  to  the 

partial  derivatives  are  second  order  and  centered  in  both  time  and  space.  The 

model  uses  a  staggered  horizontal  grid  with  variables  Ax  =  Ay  =  5  to  20  km  over 

a  2f)  y  2fi  grid-point  domain.  Velocity  components  are  defined  at  the  cross  points, 

while  all  other  variables  are  defined  at  the  dot  points  (see  Figure  1).  The  use  of 

the  staggered  grid  leads  to  a  smaller  spatial  truncation  error  of  the  pressure 

gradient  term  than  in  a  non-staggered  grid.  In  addition,  boundary  values  of  the 

3, 5 

velocity  components  do  not  appear  in  calculations  of  the  continuity  equation. 

No  map  scale  factor  has  been  incorporated  in  this  system  of  equations.  There 
arc  IT)  tci-rain -following  computational  levels  in  the  vertical.  The  distribution  of 
the  computational  levels  was  transformed  from  the  conventional  sigma  coordinate 

(I  =  (p  -  Pj.l/’f  ,  (1) 

whci'c  p  is  pressure,  p.p  is  the  pressure  at  the  top  of  the  model,  and  ‘f  =  (p^  -  Prp), 
in  which  p^  is  the  surface  pressure.  When  p,p  =  0,  a  =  p/Pg-  The  transformation 
cmeates  a  vertical  coordinate  (v)  related  to  sigma  by 

a  =  (4r  -  !/^)/3  .  (2) 

As  shown  in  Table  1,  the  levels  at  which  the  vertical  velocity  (v  =  dv/dt)  is  com¬ 
puted  are  staggered  in  the  vertical  with  respect  to  the  other  variables  in  order  to 
reduce  the  truncation  error  of  the  vertical  finite  differencing.  The  lower  bound¬ 
ary  is  defined  by  a  5'  X  5'  terrain  data  set  archived  at  the  Defense  Mapping  Agency 
and  Aeronautics  Center  (DMAAC)  in  St.  Louis,  Mo. 

The  basic  equations  are  expressed  in  flux  form  in  the  r -coordinate  system. 
The  mass -weighted  prognostic  variables  are 

U  =  TTU  , 

V  =  JTV  , 

S  =  JTs  ,  and 


5.  Anthes,  R.A,,  and  Warner,  T.  T.  (197  8)  Development  of  hydrodynamic  models 
suitable  for  air  pollution  and  other  meteorological  studies,  Mon.  Wea.  Rev. 
106:1045-1078. 
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Figure  1.  Staggered  Horizontal  Grid.  Velocity  compo¬ 
nents  are  defined  on  the  27  X  27  x-grid.  All  other  vari¬ 
ables  are  defined  on  the  26  X  26  o-grid 


where  S  is  the  entropy  (defined  below),  and  and  represent  the  mixing  ratios 
for  water  vapor  and  cloud  water  respectively.  There  are  five  predictive  equations: 

3U  3(Uu)  3(Uv)  I  3(a'Ur)  .  ,  L  dn  d{n<i>)  ,  ^ 

TT  =  ■  - 53r  ■  —57—  +  fv  +  U - p- 

'  '  (4) 

3V  3(Vu)  3(Vv)  1  3(a'Vr)  RT^ottI  3»r  d(w<t,)  ^  ^ 

■3t  "  ■  ~37J  33r  ■  —dT- - tU  -H  1  0  p— I  ^ 


3(Su)  3(Sv)  1  3(o'Sr) 


3(Wu)  _  3(Wv)  _  1  3(o'Wr)  ^  p 
3x  3y  3i/  w 


an  r  /  3U  ^  3V  \  . 

Tt  ■  -J  sr*  ay 


Table  1.  u-Levels  at  a  Gridpoint  Over  325  m  Terrain 


Level 

Pressure 

(mbar) 

(m  MSL) 

Variables  Defined 

Top 

0 

Infinite 

di//dt  =  0 

1.0 

44 

21562 

All  Other  Variables 

1.  5 

87 

di>/ dt 

2.0 

131 

14626 

All  Other  Variables 

2.5 

175 

- 

. . . etc.  . . 

3.0 

218 

11448 

3.5 

262 

- 

4.0 

305 

9251 

4.  5 

348 

- 

5.0 

391 

7509 

5.  5 

433 

- 

6.0 

475 

6074 

6.  5 

516 

- 

7.0 

556 

4871 

7.  5 

596 

- 

8.0 

635 

3847 

8.  5 

672 

- 

9.0 

709 

2971 

9.  5 

744 

* 

10.  0 

777 

2226 

10.  5 

809 

- 

11.0 

839 

1609 

11.  5 

866 

- 

12.0 

892 

1116 

12.  5 

914 

- 

13.0 

934 

740 

13.  5 

951 

- 

14.0 

965 

480 

All  Other  Variables 

14.  5 

97  5 

- 

di'/ dt 

15.0 

981 

343 

All  Other  Variables 

Bottom 

983 

325 

di/dt  =  0 

The  friction/diffusion  terms  (F^^.  F^.  F^,  and  F^)  will  be  explained  in  Section 
2.2.  Three  diagnostic  equations  provide  information  for  the  prognostic  equations: 


-  -  Cp  0(1  +  0.61  q^)  ,  and  (10a) 

5  =  7r[ln(T/P)  +  Lq^/CpT]  .  (lOo) 

The  following  definitions  apply; 
a '  =  da  /  dv  (f,  x 

=  dr/dt  T*  =  T(1  +  0.61  qy) 

P  =  (P/P^)^/S  L  =  597.  3  -  0.  566  (T  -  273.  16) 

6  =  T/P  . 

The  saturation  vapor  pressure  with  respect  to  water,  e^,  is  taken  from  Murray: 

Og  =  6.  U  exp  [  17.  27  (T  -  273.  16)/(T  -  35.  86)J  .  (11) 

In  order  to  determine  the  temperature  and  moisture  variables,  the  model 
first  computes  (saturation  temperature)  and  q^^  (saturation  mixing  ratio).  If 
the  air  is  saturated  (W  >  ^fq^g).  the  mixing  ratio  is  a  known  function  of  tempera¬ 
ture,  and  Eq.  (10b)  bec  omes  a  transcendental  equation  for  the  temperature,  T^, 
c:orresponding  to  saturation  with  respect  to  liquid  water.  If,  on  the  other  hand, 
the  air  is  unsaturated  (W  £  the  mixing  ratio  in  Eq.  (lOb)  is  replaced  by 

W/ff  and  the  temperature  is  solved  for  directly.  The  finite  difference  forms  of 

3 

Eqs.  (4)  through  (10)  are  given  by  Nickerson. 

2.2  Numerical  MetlicKls  and  I’livsical  Parainelerizalkms 

A  time  step  of  between  10  and  30  sec  may  be  used  in  this  model.  Two  differ¬ 
ent  procedures  have  been  developed  to  control  the  small  scale,  high  frequency 
waves  generated  by  the  model.  The  first  uses  the  Matsuno  Time -and -Space -Un- 
<  entered  (T.-\SU)‘  procedure  in  combination  with  a  five-point  spatial  smoother 

fi  Murray,  E.W,  (1007)  On  the  computation  of  saturation  vapor  pressure, 

I  Appl.  Meteorol.  0:203-204. 

T  Arakawa,  .A.  ,  and  Mintz,  Y.  (1974)  The  UCLA  general  circulation  model. 

Notes  from  the  workshop,  25  March  -  4  April  1974,  Department  of  Meteor¬ 
ology,  UCLA. 


applied  to  the  forecast  winds  on  r-surfaces  at  each  time  step.  The  TASU  proce¬ 
dure  requires  that  two  time  steps  (each  ninth  and  tenth)  be  broken  up  into  four 

g 

half-time  steps.  Unfortunately,  the  five-point  operator  has  been  shown  to  re¬ 
move  energy  from  the  meteorologically  significant  gravity-inertia  (GI)  wave 
modes  (as  well  as  the  non -meteorological  GI  wave  modes)  generated  during  the 
geostrophic  adjustment  process.  The  other  time -and-space  smoothing  procedure 
uses  an  Asselin  filter  in  conjunction  with  explicit,  fourth-order  horizontal  diffu¬ 
sion  terms  within  the  momentum  equations.  The  Asselin  filter  produces  some¬ 
what  smoother  time  integrations  without  the  need  for  extra  half-time  steps  as  in 
the  TASU  scheme.  However,  the  latter  requires  a  shorter  time  step  in  order  to 
maintain  computational  stability;  thus,  the  overall  computation  time  is  increased 
(by  about  20  percent  in  recent  tests).  The  diffusion  terms  and  in  Eqs.  (4) 
and  (5)  serve  two  purposes:  the  first  is  to  control  energy  growth  over  certain 
spatial  scales  in  the  model  domain;  the  second  is  to  represent  the  diffusion  proc¬ 
esses  that  actually  occur  in  the  real  atmosphere. 

2  1 1 

At  the  top  boundary,  i'  =  0.  An  upper  absorbing  layer  or  "sponge"  was  de- 

4 

veloped  and  tested.  This  layer  permits  certain  vertically-propagating  GI  waves 
(generated  by  various  local  forcing  mechanisms)  to  exit  the  model  domain. 

As  stated  in  the  introduction,  one  of  the  most  difficult  problems  with  limited- 
area  models  is  the  proper  treatment  of  the  variables  at  the  lateral  boundaries. 
This,  in  turn,  greatly  affects  the  accuracy  of  solutions  on  the  interior  of  the 
domain.  lixcept  for  outflow  winds  extrapolated  outward  from  interior  values, 

the  NOAA  model  originally  used  time -invariant  lateral  boundary  conditions. 

4  1 0 

tierently,  Nickerson  incorporated  an  Orlanski  radiation  boundary  condition 

for  the  latei-al  boundaries. 

rtie  model  treats  the  first  vertical  computational  layer  between  the  lower 
boundary  and  the  first  grid  point  above  the  surface  (approximately  18  m  thick) 
with  tlie  llusinger-Dyer  formulation^’  '  *  for  a  constant-flux  surface  layer.  This 
formulation  expresses  the  friction  veloc  ity  and  surface  fluxes  of  sensible  and 
latent  heat  as  non-dimensional  functions  (F  and  G)  of  height  and  the  Monin-Obukov 
length.  Which  particular  set  of  equations  is  used  depends  on  whether  the  surface 
layer  is  unstable,  mildly  stable,  or  stable. 

!!.  Seittcr,  K.  U.  ( 1 ‘>84)  personal  cori-espondence. 
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A  boundary  layer  is  defined  between  the  first  grid  point  above  the  surface 
(Z.,  =  18  ni)  and  the  layer  (Z.  =  1000  m),  and  eddy  exchange  coefficients  are 
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computed  in  accordance  with  the  profile  given  by  O'Brien.  Between  Z^  and 
Zg,  separate  coefficients  are  computed  for  momentum  (K^)  and  for  temperature 
and  moisture  (Krp)  by  equations  analogous  to 

K(Z)  =  K(Z^)  +  (Z  -  Z^)^/(AZ)^ 


+  (Z  -  Zp)  IK'fZjj)  +  2[K(Zg)  -  K(Zg)]/AZ 


KJZ)  =  u-/F'(Z)’ 
K„(Z)  =  u-/G'(Z) 


K^(Z)  =  0 


K.j,(Z)  =  0 


at  Z  =  Zy,  and 


for  Z  Z, 


where  u us  the  fric  tion  veloc  ity,  the  prime  denotes  the  derivative  with  respect  to 
Z,  and  AZ  =  (Z^  -  Zj^).  The  friedion  terms  in  Eqs.  (4)  through  (7)  are  written  as: 


K  =  A  ^  [  AK(u)  ^ 
u  dr  dr 


=  A  ^  AK(u) 


F  =  A  ^  (  AK(T) 

S  \  Ol/ 


F  =  A  [  AK(T)  ^ 
w  dr  dr 


ichere  A  =  -  gp/RTira'.  At  the  lower  boundary  (r  =  a  =  1), 
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AKtUj  =  m*  cos  , 

AK(u) =  irum^sin^  , 

-o  /q  P  LQ  \ 

A"™-!-  -  .  ana 

AK(T)-|^=  irQ^  ,  (14b) 

where  tan  a  =  (the  ratio  of  the  two  wind  components  at  the  first  grid  point 

above  the  surface),  and  the  surface  fluxes  of  heat  and  moisture  (Q^  and  Q^)  are 
given  by 

'Qs  =  .  and 

The  NOAA  model  employs  detailed  microphysical  parameterizations  of  warm 
rain  processes  (the  "warm"  version)  and  frozen  precipitation  processes  (the 
"cold"  version).  The  warm  version  begins  with  a  forecast  of  W  and  an  existing 
rain  water  mass  with  a  log-normal  drop  size  distribution.  Raindrops  do  not  grow 
by  diffusion,  but  may  form  and  disappear  by  separate  processes  of  evaporation, 
advection,  sedimentation,  raindrop  self-collection,  and  cloud  drop  accretion.  In 
order  to  keep  track  of  the  effect  of  these  processes,  the  model  carries  two  rain 
descriptors  at  each  grid  point;  the  rain  water  mixing  ratio  and  the  raindrop  spe¬ 
cific  number  concentration.  The  cold  version  additionally  forecasts  three  types 
of  ice  particles:  unrimed,  partially  rimed,  and  graupel.  Therefore,  six  precip¬ 
itation  variables  are  necessary  at  each  grid  point,  namely,  mixing  ratio  and 
number  concentration  for  each  of  the  three  ice  particle  types.  These  particles 
are  assumed  to  exhibit  a  Marshall-Palmer  size  distribution.  Physical  processes 
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in  the  cold  model  include  advection,  sedimentation,  spontaneous  nucleation  of 
unrimed  crystals,  accretion  of  ice  particles,  deposition/sublimation,  and  trans¬ 
fer  of  growing  ice  particles  to  the  next  larger  category.  In  both  the  warm  and 
cold  versions,  precipitation  that  falls  through  the  lowest  computational  level  ac¬ 
cumulates  on  the  ground. 

The  model  has  no  parameterizations  for  radiative  surface  heating  and  cooling, 
sub -grid  scale  convection,  or  surface  characteristics  such  as  vegetation  and  soil 
moisture.  Surface  temperature  and  moisture  are  maintained  constant  in  time. 


2.3  Model  Initialization 


A  description  of  the  initialization  procedure  of  the  model  is  given  here  to 
help  better  understand  several  problems  that  are  noted  in  Section  3. 

First,  the  model  computes  surface  pressure  at  each  grid  point  from 


”  ^o 


T  +  \®^^d”'^d 


(15) 


where  p  is  the  surface  pressure,  p  the  gridded  value  of  sea  level  pressure,  Z„ 


the  terrain  elevation,  and  T  an  assumed  mean  virtual  temperature.  A  standard 


lapse  rate  is  assumed  (y  =  0.0065°K/m),  and  g,  R  ,,  and  y,  are  gravitational 


s  -  -  — -  -  d'  'd 

acceleration,  the  gas  constant  for  dry  air,  and  the  dry  adiabatic  lapse  rate,  re¬ 
spectively.  The  winds,  temperature,  and  moisture  are  initialized  from  a  single 
radiosonde  located  at  the  grid  center.  The  model  computes  pressure  at  each  grid 
point  on  the  constant  i/ -surfaces.  The  pressure  on  each  !• -surface  is  a  constant 
multiple  of  the  underlying  surface  pressure.  Then,  temperatures  and  mixing 
ratios  are  determined  at  each  r-surface  from  a  linear-in-pressure  vertical  inter¬ 
polation  of  the  sounding  data.  An  unfortunate  consequence  of  the  single -sounding 
initialization  is  that,  on  a  given  r-surface,  horizontal  gradients  of  temperature 
and  moisture  arise  solely  from  the  existence  of  underlying  gradients  of  surface 
elevation.  The  model  insures  that  the  temperature  at  the  first  r -level  above 
ground  (level  1.5)  is  less  than  the  ground  temperature  (level  16).  Then,  the  hyp¬ 
sometric  equation  [Eq,  (10a)|  is  integrated  upward  from  the  surface  to  determine 
geopotential,  without  direct  reference  to  the  sounding.  As  a  result,  the  horizon¬ 
tal  continuity  of  geopotential  on  r -surfaces  is  also  affected  by  the  underlying 
gradients  of  surface  elevation.  Finally,  the  u-  and  v-wind  components  are  simi¬ 
larly  interpolated  from  the  sounding  to  the  r -surfaces.  All  other  model  variables 
(e.g.  ,  entropy)  are  computed  from  these  gridded  variables. 
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The  model  is  initialized  with  wind  fields  out  of  dynamic  balance  with  the 
mass  distribution.  Previous  tests  by  Nickerson  and  others  have  suggested  that 
an  elapsed  time  of  one  model-hour  is  needed  for  the  mutual  adjustment  of  these 
fields . 

2.4  Model  Configurations  Under  Development 

At  this  time,  there  are  two  different  versions  of  the  model  under  develop¬ 
ment.  The  first  is  the  AFGL/NOAA  Mesoscale  Model,  described  in  Table  2. 

The  model  is  written  in  ANSI  FORTRAN  77.  It  requires  70,  592  decimal  words 
of  central  memory,  which  economizes  on  large -core -memory  requirements  at 
the  expense  of  run  time.  In  this  version,  the  numerous  26  X  26  X  16  point  vari¬ 
able  arrays  generally  reside  in  mass  storage.  One  vertical  east-west  plane  is 
buffered  into  core  at  a  time;  hence,  there  is  never  more  than  about  15  percent  of 
the  data  in  core  at  once.  Run  time  for  a  7-h  integration  with  30-sec  time  steps 
on  the  CDC  Cyber  170/750  includes  200  min  central  processing  time  and  790  min 
input/output  time. 


Table  2.  Features  of  the  AFGL/NOAA  Model  and  the  NOAA/LAMP  Model 


Feature 

AFGL/NOAA 

NOAA/LAMP 

Microphysics 

Warm  or  Cold 

Cold 

Time  Integration 

Centered  Differences 

Centered  Differences 

Time/Space 

Smoothing 

Matsuno  TASU/5 -Point 
Smoother 

Asselin/Explicit 

Diffusion 

Upper  Boundary 
Conditions 

Constant  (Reflective) 

Sponge  or  Constant 

Lateral  Boundary 
Conditions 

Time -Invariant/Extra¬ 
polate  Outflow  Wind 

Orlanski  Radiative 

Vertical  Levels 

15 

15  (31  in  the  2D 

Version) 

Computer  Use 
Strategy 

15  Percent  Core, 

85  Percent  Mass 

Storage 

Fully  Core -Resident 

Computer 

CDC  Cyber  170/750 

Cray  1 

The  second  version  is  under  development  at  NOAA/ERL  and  in  France  by  co¬ 
operative  agreement  with  the  Laboratoire  Associe  de  Meteorologie  Physique 
(LAMP)  at  Clermont.  The  NOAA/LAMP  model  is  described  in  Table  2.  It  is 
frequently  useful  to  test  numerical  aspects  of  a  model,  such  as  vertical  resolution 
or  the  parameterization  of  physical  processes,  with  a  less  complex  and  more 


efficient  2D  analog  to  the  3D  model.  Most  of  the  above  features  are  also  present 
in  a  2D  (x-z)  version  of  the  NOAA/LAMP  model. 


3.  KF-SULTS  OF  A  BASELINE  TEST 


While  the  primary  objective  during  the  first  year  of  the  AFGL  program  was 
to  perform  a  baseline  test  of  the  NOAA  3D  model,  three  secondary  achievements 
were  also  accomplished.  First,  the  model  was  adapted  to  the  AFGL  computer 
(hence  the  name  AFGL/NOAA  model).  The  warm  version  was  selected  for  the 
baseline  test.  Second,  a  baseline  data  set  was  selected  and  prepared.  A  single 
vertical  sounding  was  used  to  initialize  the  model  variables;  however,  the  base¬ 
line  data  set  also  contains  a  dense,  3 -hourly  radiosonde  set  (27  March  1982), 
along  with  surface,  rain  gauge,  and  satellite  data.  Third,  it  was  necessary  to 

write  software  that  processed  the  model  output  into  graphical  and  statistical 
1  .“i 

form  for  later  verification  against  the  baseline  data. 

The  baseline  test  revealed  several  major  problems  with  the  model.  These 
problems  appeared  to  stem  from  the  initialization  and  the  lateral  boundary  condi¬ 
tions  in  large  part,  and  they  disrupted  the  forecast  to  the  point  where  it  was  not 
possible  to  assess  the  performance  of  other  parts  of  the  program  (e.  g. ,  boundary 
layer  formulation  and  rain  microphysics). 

The  initialization  described  in  Section  2.3  was  found  to  be  far  from  adequate 
in  three  respects.  First,  the  resulting  model  analyses  of  temperature,  moisture, 
and  winds  have  no  realistic  horizontal  gradients  on  the  computational  surfaces. 
Using  the  analyzed  temperatures  in  the  vertical  integration  of  the  hydrostatic 
equation  yields  geopotential-height  fields  for  which  the  contours  are  roughly  par¬ 
allel  to  contours  of  terrain  height.  For  example,  notice  the  discrepancy  at  500 
mbar  between  the  observed  heights  (Figure  2)  and  the  corresponding  single¬ 
sounding  initialized  heights  (Figure  3);  the  direction  of  the  gradient  is  completely 
reversed  Sei-ond,  the  temperature  analysis  at  ♦he  lowest  level  is  not  represen¬ 
tative  of  the  mean  surface  and  boundary  layer  conditions  Third,  the  period  of 
the  initial  geostrophic  adjustment  process  was  excessively  long. 

Two  approaches  are  under  development  at  NOAA/ERL  and  AFGL  to  improve 
the  initialization  scheme.  They  are:  (1)  revision  of  the  single -sounding  technique 
by  imposing  a  variational  constraint  on  the  wind  and  mass  fields,  with  sensitivity 
to  boundary  layer  structure;  and  (2)  the  use  of  more  initial  data  (e.  g. ,  multiple 
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soundings,  surface  data,  or  output  from  a  larger  model).  Some  simple  experi¬ 
ments  for  improving  the  model  initialization  in  the  lower  layers  using  surface 
data  will  be  discussed  in  Appendix  A. 

As  mentioned  above,  a  major  problem  noted  in  the  baseline  test  was  the  per¬ 
sistent  imbalance  between  wind  and  mass  fields  during  the  forecasts.  The  3-h 
forecast  of  a  500  mbar  height  field  shown  in  Figure  4  normally  implies  a  north¬ 
easterly  geostrophic  wind;  yet.  the  corresponding  forecast  wind  direction  in 
Figure  5  is  westerly  to  north-westerly.  The  forecast  shows  little  evidence  of 
geostrophic  balance  after  3  h  of  integration.  It  appears  that  the  only  practical 
solution  to  this  djmamic  imbalance  problem  is  to  begin  with  a  more  realistic,  or 
balanced,  initialization. 

Most  of  the  forecast  fields  (with  the  exception  of  wind,  which  was  specifically 
smoothed  at  each  time  step)  exhibited  significant  2Ax  noise.  These  waves  did  not 
appear  very  large  in  amplitude  (see  Figure  4).  Waves  in  the  500  mbar  tempera¬ 
ture  field  (not  illustrated)  were  about  1“C  in  amplitude.  The  wave  pattern  ap¬ 
peared  to  originate  near  the  eastern  (downwind)  boundary  and  propagate  toward 
the  center  of  the  domain  during  the  forecast.  In  addition,  a  gradually  amplifying 
disturbance  was  fixed  along  the  western  boundary.  This  evidence  strongly  sug¬ 
gests  that  these  waves  probably  arose  from  inadequate  numerical  treatment  of  the 
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Figure  5.  Three -Hour  Forecast  of  500  mbar  Wind. 
Longest  vector  represents  33  msec'^  (baseline 
model  winds  were  initially  uniform  west-southwest) 


lateral  boundaries.  Figure  6  illustrates  the  sea  level  pressure  forecast.  Sur¬ 
face  pressure  is  maintained  constant  in  time  at  the  lateral  boundaries,  but  may 
vary  in  the  interior.  Notice  the  2Ax  noise  in  this  field  at  the  west  and  east 
boundaries.  One  variable  that  reacted  strongly  to  this  was  the  vertical  velocity, 
w  (Figure  7),  which  had  an  upward  maximum  of  12  cm  sec*^  near  the  eastern 
edge  and  an  even  stronger  downward  region  on  the  western  edge.  Note  the  cor¬ 
responding  precipitation  maximum  along  the  north  and  east  edges  (Figure  8). 

The  boundary  condition  problem  is  likely  to  be  solved  only  with  a  more  com¬ 
plex  numerical  treatment.  This  was  determined  after  the  model  continued  to  ex¬ 
hibit  excessive  Gl-wave  activity  even  after  it  was  initialized  with  the  actual  AVE/ 
VAS  "truth"  data.  In  this  experiment,  the  four  model  input  fields  of  ft,  q^,  u, 
and  V  were  objectively  analyzed  to  the  model  grid  points  (see  Figure  9).  No  bal¬ 
ancing  of  the  mass  and  momentum  fields  was  performed.  The  model  was  then 
run  for  five  hours,  and  output  fields  were  compared  with  the  previously  discussed 
model  run.  While  the  result  of  this  latter  initialization  may  be  viewed  as  a  sig¬ 
nificant  improvement  from  the  single -sounding  case.  Figure  10  (500  mbar  geo¬ 
potential  height  in  m)  shows  that  large  GI  waves  continue  to  dominate  the  numer¬ 
ical  solution.  Higher  frequency  noise  also  continues  to  be  a  problem  near  the 
boundaries,  suggesting  that  the  time -invariant  boundary  conditions  may  be  a 
source  of  much  GI  energy  generation. 
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Figure  6.  Three -Hour  Forecast  of  Sea -Level  Pres¬ 
sure  (mbar).  Isobars  at  0.25  mbar  intervals 


Figure  7,  Three-Hour  Forecast  of  Vertical  Velocity 
(i  m  sec"').  Solid  lines  are  upward,  dashed  lines  are 
downwa  rd 


Figure  8.  Three -Hour  Forecast  of  Total  Precipi 
tation  (mm).  Contours  at  0.  5  mm  intervals 


H'ijfurt'  10.  Five-Hour  Forecast  of  500  mbar 
Hcif'ht  (m) 


t.  C.ONCI.LSIONS 

A  numeric  al  mesoscale  forecast  model  was  selected  for  potential  use  by  the 
Air  Force  in  an  operational  mode,  either  in  a  nested  configuration,  or  as  a  stand¬ 
alone  model  in  a  base -weather  station.  This  NOAA  model  and  many  of  its  numer¬ 
ical  and  meteorological  features  were  described  in  this  report.  Two  versions 
were  specified:  the  first  is  under  development  at  NOAA/ERL  (the  NOAA/LAMP 
Model),  and  the  other  is  under  development  at  AFGL  (the  AFGL/NOAA  Model). 

A  baseline  test  of  the  model  revealed  several  problems.  The  initialization 
Scheme  computed  initial  fields  of  mass  and  momentum  that  failed  to  mutually  ad- 
ju.st  after  several  hours  of  integration.  Work  is  under  way  at  AFGL  and  at  ERL 
to  improve  this  initialization:  subsequent  techniques  designed  to  improve  the  ini¬ 
tial  conditions  appeared  to  be  a  step  in  the  right  direction.  The  nonmeteorologi- 
cal  GI  wave  noise  evident  in  most  of  the  forecast  fields  appeared  to  result  from 
poor  treatment  of  the  variables  at  the  lateral  boundaries. 
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Appendix  A 

Alternative  Boundary  Layer  Initialization  Techniques 

Initializing  a  mesoscale  model  with  available  conventional  surface  and  upper- 
air  data  can  be  a  difficult  problem.  Sufficient  surface  observations  are  generally 
available  to  obtain  a  good  set  of  initial  surface  analyses.  However,  upper-air 
observations  are  nearly  always  far  too  sparse  to  describe  the  details  of  atmos¬ 
pheric  structure,  especially  in  the  boundary  layer.  The  AFGL  mesoscale  model, 
with  20-km  horizontal  resolution,  encompasses  a  500  500  km  area.  Within  this 

area,  there  are  in  excess  of  25  surface  stations,  but  only  one  upper-air  station. 
This  distribution  of  observations,  which  is  typical  across  the  United  States,  sug¬ 
gests  that  radiosonde  data  may  be  supplemented  with  data  extracted  from  the  sur¬ 
face  so  that  sufficient  detail  can  be  obtained  to  ac:curately  initialize  a  mesoscale 
model. 

The  current  version  of  the  AFGL  mesoscale  model  uses  a  single  sounding  to 
initialize  model  fields  at  all  levels  over  the  entire  model  domain.  This  has  been 
shown  to  be  relatively  successful  in  certain  situations  where  local  terrain  forcing 
dominates  the  atmospheric  flow.  However,  in  regions  where  terrain  forcing  is 
less,  the  initial  horizontal  fields  generated  by  the  single -sounding  model  initial¬ 
ization  tend  to  reflec  t  the  terrain  rather  than  the  actual  distribution  of  atmos¬ 
pheric  parameters.  These  single-sounding  model-initialized  fields,  therefore, 
are  physic  ally  unrealistic,  and,  consequently,  a  model  initialized  with  these 
fields  is  unlikely  to  generate  useful  forecasts. 

■As  a  first  attempt  to  improve  the  initialization  of  the  .AFGL  mesoscale  model, 
several  simple  techniques,  using  the  relatively  dense  surface  observations  in 
conjunction  with  the  <  urrent  single -sounding  model  initialization  scheme,  were 
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developed  to  generate  a  "blended"  initialization  for  the  lower  model  levels.  The 
idea  was  that  the  detail  and  realism  of  the  lower  atmospheric  fields  used  to  ini¬ 
tialize  the  model  could  be  greatly  enhanced  by  the  valuable  information  contained 
in  the  surface  analysis.  The  blended  fields  were  calculated  for  model  levels  15, 

14,  13,  12,  and  11,  approximately  corresponding  to  pressure  levels  of  998,  981, 
950,  907,  and  853  mbar,  respectively.  These  model  levels  were  chosen,  because 
analyzed  horizontal  fields  of  quantities  generally  resemble  surface  fields. 

The  fields  necessary  to  initialize  the  model  are  potential  temperature,  spe¬ 
cific  humidity,  and  the  U  and  V  components  of  the  horizontal  wind.  The  initial 
fields  generated  by  the  blending  technique  and  the  single -sounding  model  initial¬ 
ization  fields  were  then  compared  to  a  "truth  set"  of  analyzed  radiosonde  obser¬ 
vations.  The  truth  data  set  was  generated  from  a  network  of  2  1  densely  spaced 
AVE/VAS  and  National  Weather  Service  radiosonde  sites  located  in  Texas  and 
Oklahoma.  The  truth  data  set,  available  at  25-mbar  intervals,  was  interpolated 
to  the  model  r-levels.  The  model  r-level  data  was  then  analyzed  to  a  17  X  16 
(l/2  deg)  grid  using  a  2 -pass  Barnes  analysis  routine.  The  1/2  deg  data  was  then 
horizontally  interpolated  to  the  smaller  26  X  26  20-km  model  grid. 

The  truth  data  set  was  then  compared,  grid  point  by  grid  point,  to  the  single - 
sounding  model  initialization  fields.  Root-Mean-Square  Error  (RMSE)  and  Mean 
Absolute  Error  (MAE)  statistics  were  computed  for  each  parameter,  at  each  v- 
level.  These  error  statistics  were  then  used  as  the  basis  in  comparing  the  im¬ 
provement  (or  lack  of  improvement)  of  the  blending  techniques  to  the  single -sound¬ 
ing  model  initialization.  Five  tests  were  performed,  and  each  is  briefly  de¬ 
scribed  in  the  following  paragraphs. 

Test  I  was  designed  to  blend  the  surface  fields  with  the  single -sounding  model 
initialization  such  that  the  blended  fields  would  be  composed  of  100  percent  sur¬ 
face  data  at  the  model  surface  (r -level  16),  and  would  decrease  according  to  a 
mass -weighted  distribution  to  0  percent  surface  data  and  100  percent  single- 
sounding  model  initialization  at  r -level  10.  The  weights  for  each  level  are  shown 
in  Table  Al.  This  would  allow  the  surface  analysis  to  influence  the  initial  model 
fields  up  to  approximately  the  850  mbar  level. 


Table  Al.  Blended-Initialization  Weights  for  Test  1 


1-15 

c-14 

1/  - 1 3 

i'- 12 

1^-11 

0.  99 

0.  91 

0.75 

0.  54 

0.26 

0.  01 

0.09 

0.25 

0.46 

0.74 

Test  II  was  designed  to  measure  the  effect  of  increasing  the  weight  given  to 
the  surface  analysis.  The  weights  were  chosen  to  be  strictly  mass  weighted, 
assuming  a  surface  pressure  of  approximately  1000  mbar.  The  weights  for  each 
level  are  shown  in  Table  A2. 

Table  All.  Blended-Initialization  Weights  for  Test  11 


Test  III  was  designed  simply  to  show  the  effect  of  ignoring  the  single-sound¬ 
ing  model  initialization  by  applying  the  actual  surface  fields  at  the  tested  model 
levels.  Thus,  the  analyzed  surface  fields  were  given  a  weight  of  1.00,  and  the 
single -sounding  model  initialization  fields  were  given  a  weight  of  0.00,  both  at 
all  I'-levels. 

Test  IV  was  designed  to  incorporate  the  vertical  distribution  of  each  param¬ 
eter  observed  from  the  single  sounding  with  the  analyzed  surface  fields.  The 
lapse  rate  of  each  parameter  was  calculated  from  the  single  sounding.  The  lapse 
rate  of  each  parameter  was  then  applied  to  the  analyzed  surface  fields  to  calcu¬ 
late  blended  fields  for  each  parameter  at  each  t/ -level. 

Test  V  was  an  attempt  to  improve  upon  the  other  tests,  all  of  which  failed  to 
generate  accurate  blended  wind  fields,  and,  thus,  was  applied  to  U  and  V  fields 
only.  An  Ekman  layer  was  assumed,  and  the  horizontal  wind  components  were 
calculated  for  each  r -level  and  blended  with  the  single -sounding  model  initializa¬ 
tion  horizontal  wind  components  with  the  same  weights  as  in  Test  I. 

Table  A3  indicates  that  Test  I  and  Test  IV  provide  significant  improvement 
over  the  single -sounding  model  initialization  fields  for  potential  temperature  and 
mixing  ratio  at  all  levels.  Test  IV  shows  the  most  improvement  overall,  exhibit¬ 
ing  the  smallest  errors  at  the  higher  levels  of  all  tests.  It  shows  improvements 
for  potential  temperature  of  77,  77,  71,  55,  and  53  percent  over  the  single¬ 
sounding  model  initialization  for  v -levels  15,  14,  13,  12,  and  11,  respectively, 
and  improvements  of  69,  65,  41,  39,  and  57  percent  at  the  same  i/ -levels  for 
mixing  ratio.  Tests  II  and  III  also  show  considerable  improvement  over  the 
single -sounding  model  initialization  for  potential  temperature  and  mixing  ratio  at 
r-levels  15,  14,  and  13;  however,  at  i/ -levels  12  and  11,  these  methods  do  poorly 
due  to  the  large  weights  given  to  the  surface  analysis. 


Table  A3.  Root-Mean-Square  Errors  of  Blended  Sounding  Technique 


While  these  blending  techniques  show  great  potential  for  improvement  of 
model  initial  analyses  over  the  single -sounding  model  initialization  technique  for 
potential  temperature  and  mixing  ratio,  the  results  for  the  horizontal  wind  com¬ 
ponents  are  quite  poor.  Table  A3  shows  some  improvement  over  the  single - 
sounding  model  initialization  in  the  lowest  layers;  however,  the  errors  for  all 
techniques  tend  to  increase  with  height  and  are  generally  worse  than  the  single¬ 
sounding  model  initialization.  Each  technique  crudely  allows  a  turning  of  the 
wind  with  height  by  including  both  a  weighted  surface  wind  and  a  weighted  single¬ 
sounding  model  initialized  wind  in  the  U  and  V  calculations,  or  by  the  lapse  rate 
or  Ekman  spiral  techniques.  They  all  ignore  stability,  terrain,  and  local  effects, 
however,  and  Tests  I -IV  assume  that  the  one  sounding  used  in  the  single -sounding 
model  initialization  is  representative  of  the  entire  model  domain. 

The  results  of  these  tests  are  based  on  only  one  case,  and,  although  they 
show  the  potential  for  significant  improvement  over  the  single -sounding  model 
initialization  for  potential  temperature  and  mixing  ratio,  they  are  not  conclusive. 
Future  tests  will  be  conducted  using  analyses  of  three  or  more  radiosonde  obser¬ 
vations  from  the  area  surrounding  the  model  domain.  The  radiosonde  observa¬ 
tions  will  be  from  the  conventional  network  of  radiosonde  sites,  and,  thus,  typi¬ 
cal  of  day-to-day  conditions.  Interpolating  these  analyses  to  the  smaller  model 
grid  will  yield  a  more  realistic  representation  of  the  atmosphere  than  the  single - 
sounding  model  initializations,  with  reasonable  gradients.  These  more  realistic 
analyses  will  then  be  blended  with  the  surface  analyses  as  in  Tests  I -IV.  More 
sophisticated  blending  techniques  must  be  developed  to  improve  the  U  and  V 
initialization  fields. 


